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V\  ^  ABSTRACT 

■'  A  5  i  '-  ' 

We  use  the  method  of  probability  weighted  moments  to  derive  estimators  of 
-the  parameters  and  quantiles  of  the  generalized  extreme-value  distribution. 

■  i  ■  t 

We  investigate  the  properties  of  these  estimators  in  large  samples,  via 
asymptotic  theory,  and  in  small  and  moderate  samples,  via  computer  simulation. 
Probability  weighted  moment  estimators  have  low  variance  and  no  severe  bias, 
and  compare  favourably  with  estimators  obtained  by  the  methods  of  maximum 
likelihood  or  sextiles.  The  method  of  probability  weighted  moments  also 
yields  a  convenient  and  powerful  test  of  whether  an  extreme-value  distribution 
is  of  Fisher-Tippett  type  I,  II  or  III. 
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SIGNIFICANCE  AND  EXPLANATION 


Many  problems  in  hydrology  and  civil  engineering  are  related  to  the 
properties  of  extreme  events:  maximum  levels  or  flow  rates  of  rivers, 
greatest  wave  heights  or  sea  levels,  or  maximum  wind  loads  on  buildings.  The 
magnitudes  of  extreme  events  are  random  quantities  whose  distribution  is  often 
well  described  by  the  generalized  extreme-value  distribution. 

The  method  of  probability  weighted  moments  is  a  general  procedure  for 
estimating  parameters  and  quantiles  of  probability  distributions.  It  is  here 
applied  to  the  generalized  extreme-value  distribution.  We  derive  the  large- 
sample  distribution  of  estimators  obtained  by  the  method  of  probability 
weighted  moments,  and  compare  the  performance  of  these  estimators  in  small 
samples  with  the  performance  of  the  currently  widely-used  maximum-likelihood 
and  sextile  estimators.  The  probability  weighted  moment  estimators  generally 
have  the  best  performance. 

We  also  consider  the  problem  of  testing  hypotheses  about  the  shape  para¬ 
meter  of  the  generalized  extreme-value  distribution  and  show  that  the  method 


of  probability  weighted  moments  yields  a  convenient  and  powerful  testing 


procedure. 
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ESTIMATION  OF  THE  GENERALIZED  EXTREME-VALUE  DISTRIBUTION 
BY  THE  METHOD  OF  PROBABILITY  WEIGHTED  MOMENTS 

J.  R.  M.  Hosking1,  J.  R.  Wallis2  and  E.  F.  Wood3 
Introduction 

The  generalized  extreme-value  distribution  of  Jenkinson  (19S5)  is  widely  used  for 
modelling  extremes  of  natural  phenomena,  and  is  of  considerable  hydrological  importance, 
since  it  is  recommended  by  the  Flood  Studies  Report  (NERC  197S)  for  modelling  the 
distribution  of  annual  maxima  of  daily  streamflows  of  British  rivers.  Currently  favoured 
methods  of  estimation  of  the  parameters  and  quantiles  of  the  distribution  are  Jenkinson' s 
(1969)  method  of  sextiles,  and  the  method  of  maximum  likelihood  (Jenkinson  1969;  Prescott 
and  Walden  1980,  1983),  Neither  method  is  completely  satisfactory:  the  justification  of 
the  maximum-likelihood  approach  is  based  on  large-sample  theory,  and  there  has  been  little 
assessment  of  the  performance  of  the  method  when  applied  to  small  or  moderate  samples; 
while  the  sextile  method  involves  an  inherent  arbitrariness  (why  sextiles  rather  than,  say, 
quartiles  or  octiles?),  requires  interpolation  in  a  table  of  values  of  a  function  in  order 
to  estimate  the  shape  parameter  of  the  distribution,  and  has  statistical  properties  which 
are  not  known  even  for  large  samples. 

Probability  weighted  moments,  a  generalization  of  the  usual  moments  of  a  probability 
distribution,  were  introduced  by  Greenwood  et  al  (1979).  There  are  several  distributions, 
for  example  the  Gumbel,  logistic  and  Weibull,  whose  parameters  can  be  conveniently  esti¬ 
mated  from  their  probability  weighted  moments.  The  Gumbel  distribution,  being  a  special 
case  of  the  generalized  extreme-value  distribution,  is  of  particular  interest.  Landwehr  et 
al  (1979)  investigated  the  small-sample  properties  of  probability  weighted  moment 
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estimators  (PWM  estimators)  of  parameters  and  quantiles  for  the  Gumbei  distribution  and 
found  them  superior  in  many  respects  to  the  conventional  moment  and  maximum- likelihood 
estimators-  Die  estimators  used  by  Iandwehr  et  al  (1979)  are  identical  to  Downton's 
(1966b)  "Linear  estimates  with  linear  coefficients",  and  thus  share  the  asymptotic 
properties  of  the  latter:  in  particular,  the  asymptotic  efficiencies  of  the  PWM  estimators 
of  the  Gumbei  scale  and  location  parameters  are  0.756  and  0.996  respectively. 

In  the  present  work  we  summarize  some  theory  for  probability  weighted  moments  and  show 
that  they  can  be  used  to  obtain  estimates  of  parameters  and  quantiles  of  the  generalized 
extreme-value  distribution.  We  derive  the  asymptotic  distributions  of  these  estimators, 
and  compare,  via  computer  simulation,  the  small-sample  properties  of  the  probability 
weighted  moment,  sextile,  and  maximum-likelihood  estimators.  The  method  of  probability 
weighted  moments  outperforms  the  other  methods  in  many  cases  and  will  usually  be  preferred 
to  them.  We  also  derive  from  the  PWM  estimator  of  the  shape  parameter  of  the  generalized 
extreme-value  distribution  a  test  of  whether  this  shape  parameter  is  zero,  and  assess  the 
performance  of  this  test  by  computer  simulation. 

Probability  weighted  moments 

The  probability  weighted  moments  of  a  random  variable  X  with  distribution  function 
F(x)  ■  P(X  <  x)  are  the  quantities 

M_  -  B(XP{F<X)}r{l  -  F(X)}8)  ,  (1) 

P»*  I* 

where  p,  r  and  s  are  real  numbers.  Probability  weighted  moments  are  likely  to  be  most 
useful  when  the  inverse  distribution  function  x(F)  can  be  written  in  closed  form,  for 
then  we  may  write 

Mp,r,s  “  J^x(F)}PFr(1-F)SdF  ,  (?) 

and  this  is  often  the  most  convenient  way  of  evaluating  these  moments.  The  quantities 
M  .  „,  p  -  1,2,...,  are  the  usual  noncentral  moments  of  X.  The  moments  M.  may, 

however,  be  preferable  for  estimating  the  parameters  of  the  distribution  of  x,  since  the 
occurrence  of  only  the  first  power  of  X  in  the  expression  for  M.  ,  .  means  that  the 
relationship  between  parameters  and  moments  often  takes  a  simpler  form  in  this  case  than 
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when  using  the  conventional  moments.  When  r  and  s  are  integers,  F^l-P)8  may  be 

expressed  as  a  linear  combination  either  of  powers  of  F  or  powers  of  (1-P),  so  it  is 

natural  to  summarize  a  distribution  either  by  the  moments  r  ”  0,1,2,...,  or  by 

M1  0  s'  3  “  0,1,2,...  .  Greenwood  et  al  (1979)  generally  favour  the  latter  approach  but 

here  we  will  consider  the  moments  8  ”  M  .  *  E(X{F(X)}r],  r  -  0,1,2,...  . 

r  i ,  r ,  u 

Given  a  random  sample  of  size  n  from  the  distribution  F,  estimation  of  8^  is 

most  conveniently  based  on  the  ordered  sample  Xj  <  x ^  <•••<  x^.  The  statistic 

H  „  -1  \  (j- 1 ) ( j-2) .  . . ( j-r )  ... 

r  (n-1)(n-2)...(n-r)  j 

is  an  unbiased  estimator  of  6  (Landwehr  et  al  1979).  Instead  one  may  estimate  flr  by 

*  -1  n  r 

8  [pj  *  n  l  <<> 

3  J=1  3  3 

where  p^  is  a  plotting  position,  i.e.  a  distribution-free  estimate  of  F(x^).  Reasonable 
choices  of  p^,  such  as  p^  -  (j-a)/n,  0  <  a  <  1,  or  Pj  ”  (J-a)/(n+1-2a),  -  ^  <  a  <  -j, 
yield  estimators  ®rtPjl  which  are  asymptotically  eguivalent  to  bf  and  therefore  are 
consistent  estimators  of  fi^. 

The  estimators  br  are  closely  related  to  U-statistlcs  (Hoeffding  1948),  which  are 
averages  of  statistics  calculated  from  all  subsamples  of  size  j  <  n  of  a  given  sample  of 
size  n.  In  particular  bQ  »  n"'ix^  is  a  trivial  example  of  a  U-statlstic,  and  is  a 
natural  estimator  of  location  of  a  distribution; 

2b 


>1  -  bo  “  \  °2  “  1  (2)'1  ^ 

i>j  3 


(5) 


is  a  U-statistic  for  estimating  the  scale  of  a  distribution  -  the  statistic  Uj,  sometimes 
known  as  Gini’s  mean  difference,  has  a  history  going  back  at  least  as  far  as  von  Andrae 
(1872),  and  -j  /w  U ^  is  a  98%  efficient  estimator  of  the  scale  parameter  of  a  Normal 


distribution  (Downton  1966a);  and 


6b2  -  6b1  +  bQ 


1V3  "3  (3)_1  ^  (Xi"2VXk) 

3333  i>j>k  3 


(6) 


is  a  U-statlstic  for  estimating  skewness  which  has  been  used  as  the  basis  of  a  test  for 
Normality  by  Locke  and  Spurrier  (1976).  U-statistics  are  widely  used  in  nonparametric 
statistics  (see,  for  example,  Fraser,  1957,  Chap.  4;  Randles  and  Wolfe,  1979,  Chap.  3)  and 
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their  desirable  properties  of  robustness  to  outliers  in  the  sample,  high  efficiency  and 


asymptotic  normality  may  be  expected  to  extend  to  the  probability  weighted  moment  esti¬ 
mators  hr  and  other  quantities  calculated  from  them. 


PWM  estimators  for  the  generalized  extreme-value  distribution 
The  generalized  extreme-value  (GEV)  distribution,  introduced  by  Jenkinson  (1955), 
combines  into  a  single  form  the  three  possible  types  of  limiting  distribution  for  extreme 
values,  as  derived  by  Fisher  and  Tippett  (1928).  The  distribution  function  is 

(expt-  {1  -  k(x-E)/a}Vk],  k  ft  0  , 

(7) 

exp(-  exp{-  (x-E)/a}]  ,  k  -  0  , 

with  x  bounded  by  5  +  a/k  from  above  if  k  >  0  and  from  below  if  k  <  0.  Here  E  and 

a  are  location  and  scale  parameters  respectively,  while  the  shape  parameter  k  determines 

which  extreme-value  distribution  is  represented:  Fisher-Tippett  types  I,  II  and  III 

correspond  to  k  “  0,  k  <  0  and  k  >  0  respectively.  When  k  *  0  the  GEV  distribution 

reduces  to  the  Gumbel  distribution.  The  inverse  distribution  function  is 

f  C  +  a{l  -  (-log  F)k}/k,  y  ft  0  , 

x(F)  -  <  (8) 

E  -  a  log! -log  F) ,  k  «=  0 

The  probability  weighted  moments  of  the  GEV  distribution  for  k  ?  0  are  given  by 

Br  -  ( r+ 1 ) ” 1  [ £  +  a{  1  -  (r+1)"kr(1+k))Al  ,  k  >  -1  (9) 

(for  proof  see  Appendix  1).  When  k  <  -1,  (the  mean  of  the  distribution)  and  the  rest 

of  the  6^  do  not  exist.  From  (9)  we  have 


80  -  E  *  a(l  -  r<l+k)}/k  , 
2Bi  -  p0  ■  armkxi  -  2-k)/k 


3B2'B0  1-J-k 

2BrB0 "  i-2-k 


(10) 

(11) 

(12) 


The  PWM  estimators  E 
and  k  when  the  6 

r 


,  a,  k  of  the  parameters  are  the  solutions  of  (10) 
are  replaced  by  their  estimators  br  or  [p^ )  . 


-  (12)  for  E,  o 
To  obtain  k  we 
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must  solve  the  equation 


The  exact  solution  requires  iterative  methods,  but  because  the  function  (1-3  ) / ( 1 —2  ) 

Is  almost  linear  over  the  range  of  values  of  k,  -  <  k  <  ■j,  which  is  usually  encountered 
in  practice,  low-order  polynomial  approximations  for  k  are  very  accurate.  We  propose  the 
approximate  estimator 


7.8590c  +  2.9554c  ,  c 


2b  -b 
1  0 

3b,  -b„ 
2  0 


the  error  in  k  due  to  using  (14)  rather  than  (13)  is  less  than  0.0009  throughout  the 
1*1  * 

range  -  —  <  k  <  — .  Given  k,  the  scale  and  location  parameters  can  be  estimated 
successively  as 

,  (2b  -b  )k  ,  . 

a,  »  - Si - —  ,  5  -  bQ  +  a{T(  1+k)  -  l}/k  .  (15) 

r(1+k)(1-2-k) 

A 

Equations  (13)  and  (15),  or  their  equivalent  forms  with  hf  replaced  by  8rtp^)f  define 
the  PWM  estimators  of  the  parameters  of  the  GEV  distribution.  Given  the  estimated  para¬ 
meters,  the  quantiles  of  the  distribution  are  estimated  using  the  inverse  distribution 
function  (8). 

When  calculated  using  br  as  the  estimator  of  the  PWM  estimates  of  the  GEV  dis¬ 

tribution  satisfy  a  feasibility  criterion,  namely  that  k  >  -1  and  a  >  0  almost  surely 
(for  proof  see  Appendix  2).  This  is  clearly  a  desirable  property,  since  one  would  like 
estimates  calculated  using  a  set  of  sample  moments  to  yield  an  estimated  distribution  for 
which  the  corresponding  population  moments  exist.  We  have  not  been  able  to  prove  that  this 
feasibility  criterion  is  satisfied  when  plotting-position  estimators  8r  [p^ )  are  used,  but 
no  examples  of  the  criterion  not  being  satisfied  have  been  discovered  in  practice. 
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Asymptotic  distribution  of  HM  estimators 


When  modelling  the  properties  of  extremes  of  hydroloqicel  variables  it  rarely  occurs 
that  the  available  data  set  is  large  enough  to  ensure  that  asymptotic  large-sample  theory 
may  be  directly  applied  to  the  problem*  It  is  nonetheless  valuable  to  investigate  the 
asymptotic  properties  of  a  new  statistical  technique,  for  two  main  reasons*  First*  one  may 
seek  to  establish  the  integrity  of  the  technique,  in  the  sense  that  when  a  large  sample  is 
available,  the  new  method  should  not  be  grossly  inefficient  compared  to  an  established, 
asymptotically  optimal  method  such  as  maximum  likelihood.  Second,  asymptotic  theory  may 
provide  an  adequate  approximation  to  scorn  aspect  of  the  distribution  of  a  statistic  even  in 
quite  small  samples.  In  the  present  case  we  shall  see  that  the  variance  of  PWM  estimators 
of  parameters  and  quantiles  of  the  GEV  distribution  is  well  approximated  by  asymptotic 
theory  for  sample  sices  of  50  or  larger. 

We  consider  first  the  asymptotic  distribution  of  the  hr.  From  (3),  br  is  a  linear 
combination  of  the  order  statistics  x.,,...,xn,  and  the  results  of  Chernoff  at  al  (1967) 
may  be  used  to  prove  that  the  vector  b  -  (bQ  b,  bj)T  has  asymptotically  a  multivariate 

T 

Normal  distribution  with  mean  6  »  ( PQ  6^  P  2>  and  covariance  matrix  n  V.  The  elements 
of  v  and  details  of  the  proof  are  given  in  Appendix  3. 

The  asymptotic  distribution  of  the  FWM  estimators  of  the  GEV  parameters  follows  fro* 

A  AAA  ip 

the  preceding  result.  let  A  -  (£  a  k)  ,  •  •  IE  a  kl  ,  and  write  (13)  and  (15)  as  the 

A 

vector  equation  8  <•  f(b).  Define  the  3x3  matrix  G  «  (g^)  bY  9jj  “  Sf^/Sb^.  Then 
asymptotically  8  has  a  multivariate  Normal  distribution  with  mean  f(P)  -  9  and 
covariance  matrix  n-1GVGT  (Rao,  1973,  p.  388).  The  covariance  matrix  has  the  form 


f  2  2  n 

a  w„  a  w12  aw13 

-1 _ T  -12  2 

n  GVG  -  n  a  w12  a  *22  awJ3 


I**  13  ~23  w33 


The  Wj^  are  functions  of  k  and  have  a  complicated  algebraic  form,  but  they  can  be 
evaluated  numerically  and  are  given  in  Table  1  for  several  values  of  k.  As  k  approaches 
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-  — ,  the  variance  of  the  GEV  distribution  becomes  infinite  and  the  variances  of  the  br 
and  of  the  FWM  parameter  estimators  are  no  longer  of  order  n”^  asymptotically . 

Table  1.  Elements  of  the  asymptotic  covariance  matrix  of  the  WM 
estimators  of  the  parameters  of  the  GEV  distribution. 


k 

W11 

w12 

w13 

w22 

w23 

W33 

-0.4 

1.6637 

1.3355 

1.  1405 

1.8461 

1.  1628 

2.9092 

-0.3 

1.4153 

0.8912 

0.5640 

1.2574 

0.4442 

1.4090 

-0.2 

1.3322 

0.6727 

0. 3926 

1.0013 

0.2697 

0.9139 

-0.1 

1.2915 

0.5104 

0.3245 

0.8440 

0.2240 

0.6815 

0.0 

1.2687 

0.3705 

0.2995 

0.7395 

0.2249 

0.5635 

0.1 

1.2551 

0.2411 

0.2966 

0.6708 

0.2447 

0.5103 

0.2 

1.2474 

0.1177 

0.3081 

0.6330 

0.2728 

0.5021 

0.3 

1.2438 

-0.0023 

0.3297 

0.6223 

0.3033 

0.5294 

0.4 

1.2433 

-0. 1205 

0.3592 

0.6368 

0.3329 

0.5880 

The  asymptotic  biases  of  the  estimators  are  of  order  n  1  and  can  be  evaluated  by 
methods  similar  to  those  of  Rao  (1973,  p.  388).  The  biases,  graphed  in  Figure  1,  are 
negligible  in  large  samples  provided  that  k  >  -0.4. 

The  asymptotic  variances  of  the  estimators  are  graphed  in  Figure  2,  and  their 

asymptotic  efficiencies  in  Figure  3.  Asymptotic  efficiency  is  defined  as  the  ratio 

*  ~ 
eff(6  )  -  lim  (var  O^var  0 ^ ) 
n+*»  1 

for  each  element  0^  of  the  parameter  vector  8,  where  9^  is  the  maximum-likelihood 
estimator  of  8^  Overall  efficiency  is  the  ratio  of  the  determinants  of  the  asymptotic 
covariance  matrices  of  9  and  9.  Hie  overall  efficiency  of  the  PWM  estimators  tends  to 
zero  at  k  »  +0.5  but  for  values  of  k  not  too  far  from  zero  the  PWM  method  is  reasonably 
efficient.  Within  the  range  -0.2  <  k  <  0.2,  which  is  valid  for  many  hydrological  data 
sets,  each  PWM  parameter  estimator  has  efficiency  of  over  0.7. 


of  the  GEV  distribution 
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Corresponding  results  may  be  obtained  for  PWM  estimators  of  quantiles  of  the  GEV 
distribution.  These  are  not  presented  in  full,  but  Tables  2  and  3  give  results  for  various 
quantiles  when  k  -  -0.2,  and  for  various  values  of  X  at  the  F  -  0.98  quantile.  The 
tables  illustrate  the  main  characteristics  of  WM  quantile  estimators,  which  are:  high 
positive  bias  in  extreme  upper  tail,  arising  from  positive  bias  in  k:  high  variance  in 
upper  tail  when  k  <  0»  fair  or  high  efficiency  except  when  k  is  close  to  ±0.5. 


Table  2.  Asymptotic  bias,  variance  and  efficiency  of  PWM  estimators 
of  GEV  quantiles.  Parameters  {  •  0,  a  *  I,  k  *  -0.2. 


F 

X  (  F) 

n  x 

bias 

n  x 

variance 

efficiency 

0.001 

-1.60 

-1.2 

3.78 

0.60 

0.01 

-1.32 

-0.2 

2.06 

0.66 

0.1 

-0.77 

0.8 

0.86 

0.92 

0.2 

-0.45 

1.0 

0.88 

1.00 

0.5 

0.38 

0.6 

1.92 

0.93 

0.8 

1.75 

-1.3 

6.10 

0.98 

0.9 

2.84 

-3.1 

16.1 

0.99 

0.98 

5.91 

-4.4 

147 

0.89 

0.99 

7.55 

-1.6 

336 

0.86 

0.998 

12.33 

23.9 

1760 

0.81 

0.999 

14.90 

49.2 

3310 

0.80 
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Table  3.  Asymptotic  bias,  variance  and  efficiency  of  the  FWM  estimator 


of  the 
K  «  0, 

F  »  0.98 
a  -  1. 

quantile  of  the  GEV 

distribution 

n  x 

n  x 

k 

x(  F) 

bias 

variance 

efficiency 

-0.4 

9.41 

-64.8 

1170 

0.49 

-0.3 

7.41 

-15.3 

369 

0.75 

-0.2 

5.91 

-4.4 

147 

0.89 

-0.1 

4.77 

-1.1 

65.8 

0.96 

0.0 

3.90 

-0.1 

29.6 

0.95 

0.1 

3.23 

0.3 

14.7 

0.88 

0.2 

2.71 

0.5 

7.53 

0.75 

0.3 

2.30 

0.5 

4.04 

0.56 

0.4 

1.98 

0.6 

2.28 

0.36 

The  results  of  this  section  were  derived  for  MM  estimators  which  use  br  to  estimate 

A 

8f.  If  the  plotting-position  estimates  tp^]  are  used  instead,  the  asymptotic  variances 
and  efficiencies  remain  unchanged,  but  the  asymptotic  biases  are  different  and  cannot  be 
easily  calculated,  being  affected  by  the  biases  in  the  tp^  1  themselves. 

asall-eample  properties  of  estimates  of  the  GEV  distribution 
A  computer  simulation  experiment  was  run  to  compare  three  methods  of  estimation  of  the 
parameters  and  quantiles  of  the  GEV  distribution.  Simulations  were  perforated  for  sample 
sixes  n  «  15,  25,  50,  100  with  the  shape  parameter  of  the  distribution  taking  values  k  » 
-0.4,  -0.2,  0.0,  0.2,  0.4.  All  the  methods  of  estimation  are  Invariant  under  linear  trans- 
formatione  of  the  data,  so  without  loss  of  generality  the  location  and  scale  parameters 
C  ■  0  and  a  ■  1  were  used  throughout.  For  each  combination  of  values  n  and  k,  1000 
random  samples  were  generated  from  the  GEV  distribution,  and  for  each  sample  the  parameters 
5,  a  and  k,  and  the  quantiles  x(F),  where  F  -  0.001,  0.01,  0.1,  0.2,  0.5,  0.8,  0.9, 
0.98,  0.99,  0.998,  0.999,  were  estimated  by  each  of  three  methods!  (1)  the  method  of 
probability  weighted  moments  (PWM),  described  above!  (2)  the  method  of  maximum  likelihood 
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(ML),  using  Newton-Raphson  iteration  to  maximize  the  likelihood  function,  as  recommended  by 
Prescott  and  Walden  (1983)  and  implemented  by  Hoskinq  (1984b);  and  (3)  Jenkinson's  (1969) 
method  of  sextiles  (JS).  The  PWM  method  requires  a  choice  of  a  suitable  estimator  of  R^. 
Several  possibilities  were  investigated,  including  the  unbiased  estimator  br  and  a  number 
of  plotting-position  estimators  B^tp^].  Hie  best  overall  results  were  given  by  the 
estimator  B  (p . ]  with  p^  *  (j-0.35)/n,  and  the  simulation  results  presented  below  for 
the  PWM  method  refer  to  this  version  of  the  PWM  estimators. 

For  some  simulated  samples,  the  maximum-likelihood  and  sextile  estimates  could  not  be 
found.  The  cause  of  this  problem  for  the  maximum-likelihood  method  was  nonconvergence  of 
the  Newton-Raphson  iteration,  usually  due  to  an  extreme  outlier  in  the  sample;  for  the 
sextile  method,  a  ratio  of  sextile  means  used  to  estimate  the  shape  parameter  of  the  GEV 
distribution  sometimes  lay  outside  the  range  of  a  table  of  values  in  which  it  was  to  be 
interpolated.  Such  cases  were  omitted  from  the  simulations.  No  such  problems  were 
encountered  with  the  PWM  estimators,  which  could  always  be  calculated. 

The  simulation  results  for  estimation  of  the  parameters  of  the  GEV  distribution  are 
summarized  in  Table  4  and  5.  Results  for  the  estimator  of  k  are  of  the  greatest 
importance,  since  this  parameter  determines  the  overall  shape  of  the  GEV  distribution  and 
the  rate  of  increase  of  the  upper  quantiles  x(F)  as  F  approaches  1.  The  PWM  estimator 
has  the  lowest  standard  deviation  of  the  three  methods,  except  in  the  case  k  -  0.4,  and 
its  advantage  is  particularly  marked  in  small  samples,  n  «  15  and  n  =  25.  The  PWM  esti¬ 
mator  is  more  biased  than  the  maximum-likelihood  estimator  but  its  bias  is  small  near  the 
important  value  k  -  0.  The  sextile  estimator  of  k  has  a  large  positive  bias  in  small 
samples  when  k  <  0  and  its  standard  deviation  is  generally  larger  than  that  of  the  PWM 
estimator. 
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Table  4.  Bias  of  aatiaatora  of  GEV  parameter a 
Bias(C)  Bias (a)  Blaa(k) 


n 

ki  -0.4 

Method 

-0.2 

0.0 

0.2 

0.4 

-0.4 

-0.2 

0.0 

0.2 

0.4 

-0.4 

-0.2 

0.0 

0.2 

0.4 

MM 

.06 

.02 

-.01 

-.04 

-.  06 

.00 

-.  06 

-.09 

-.10 

-.10 

.11 

.04 

.02 

-.07 

-.11 

15 

ML 

.00 

.01 

.01 

.01 

.00 

-.08 

-.07 

-.05 

-.04 

-.  03 

-.05 

-.  02 

.02 

.03 

.04 

JS 

.29 

.25 

.22 

.17 

.10 

.23 

.12 

.06 

.00 

-.05 

.18 

.13 

.11 

.07 

-.02 

MM 

.06 

.03 

.01 

.00 

-.02 

.01 

-.04 

-.06 

-.07 

-.07 

.07 

.02 

-.02 

-.05 

-.08 

25 

ML 

.02 

.02 

.03 

.04 

.04 

-.04 

-.04 

-.03 

-.  03 

-.  02 

-.02 

-.01 

.01 

.03 

.05 

JS 

.12 

.12 

.10 

.09 

.07 

.13 

.03 

-.01 

-.04 

-.08 

.10 

.06 

.04 

.02 

-.01 

MM 

.04 

.02 

.01 

.00 

-.01 

.01 

-.02 

-.03 

-.03 

-.03 

.05 

.02 

.00 

-.02 

-.04 

SO 

ML 

.01 

.01 

.01 

.02 

.02 

-.02 

-.02 

-.  02 

-.01 

-.01 

-.01 

.00 

.01 

.02 

.03 

JS 

.05 

.03 

.02 

.01 

.01 

.04 

.01 

-.01 

-.02 

-.02 

.07 

.04 

.02 

.01 

.00 

MM 

.01 

.00 

.00 

.00 

.00 

.00 

.00 

.00 

.00 

.00 

.03 

.01 

.00 

.00 

.00 

100 

ML 

.00 

.00 

.00 

.00 

.01 

.01 

.01 

.01 

.01 

.00 

.00 

.00 

.00 

.01 

.02 

JS 

.04 

.03 

.03 

.03 

.03 

.10 

.04 

.03 

.01 

.00 

.05 

.02 

.02 

.02 

.02 

Tabla  S.  Standard  deviation  of  aatiaatora  of  GEV  parameter* 

*  *  * 
St. dev. (Cl  St.dev(a)  St.dev(k) 


n 

k«  -0.4 
Method 

-0.2 

0.0 

0.2 

0.4 

-0.4 

-0.2 

0.0 

0.2 

0.4 

-0.4 

-0.2 

0.0 

0.2 

0.4 

MM 

.30 

.29 

.28 

.28 

.28 

.32 

.25 

.20 

.18 

.18 

.21 

.20 

.20 

.20 

.21 

15 

ML 

.33 

.32 

.31 

.30 

.29 

.27 

.25 

.23 

.22 

.22 

.38 

.36 

.34 

.31 

.29 

JS 

.37 

.33 

.30 

.28 

.25 

.43 

.31 

.26 

.21 

.18 

.25 

.25 

.24 

.21 

.17 

MM 

.25 

.24 

.23 

.23 

.23 

.25 

.20 

.17 

.15 

.15 

.17 

.16 

.14 

.14 

.15 

25 

ML 

.24 

.24 

.24 

.24 

.23 

.21 

.19 

.17 

.16 

.17 

.24 

.21 

.19 

.17 

.17 

JS 

.28 

.25 

.24 

.23 

.22 

.32 

.21 

.18 

.15 

.14 

.17 

.18 

.17 

.16 

.14 

MM 

.17 

.17 

.16 

.16 

.16 

.17 

.14 

.12 

.11 

.11 

.14 

.12 

.11 

.10 

.11 

50 

ML 

.17 

.17 

.17 

.16 

.16 

.15 

.13 

.12 

.11 

.11 

.15 

.13 

.12 

.11 

.11 

JS 

.18 

.17 

.16 

.16 

.16 

.16 

.14 

.12 

.11 

.11 

.15 

.13 

.11 

.11 

.11 

MM 

.12 

.12 

.12 

.11 

.11 

.12 

.10 

.09 

.08 

.08 

.11 

.09 

.07 

.07 

.07 

00 

ML 

.12 

.12 

.12 

.12 

.11 

.10 

.09 

.08 

.08 

.08 

.10 

.09 

.08 

.07 

.06 

JS 

.12 

.12 

.12 

.12 

.11 

.16 

.11 

.09 

.08 

.08 

.09 

.09 

.08 

.08 

.08 
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Similar  results  can  b«  seen  for  eatimatora  of  £  and  a.  In  general .  MM  eatimatora 
have  amalleat  standard  deviation,  particularly  for  n  -  15  and  n  “  25,  and  their  biaa  ia 
not  large.  The  standard  deviations  of  the  MM  eatimatora  for  n  >  50  are  well 
approximated  by  their  large-sample  values  given  by  (16)  and  Table  1.  Maximum-likelihood 
estimators  are  the  least  biased  but  are  more  variable  than  MM  estimators  in  small 
samples.  Even  at  sample  sire  100,  the  asymptotic  inefficiency  of  the  PWt  method  compared 
to  maximum  likelihood  is  not  apparent  in  the  simulation  results.  Sextile  estimators  in 
general  have  larger  standard  deviations  than  PW  estimators  and  have  some  large  biases  in 
small  samples  when  k  <  0. 

The  statistical  properties  of  estimators  of  quantiles  of  the  GEV  distribution  were 
evaluated  for  many  combinations  of  quantiles  and  values  of  the  shape  parameter  k,  and 
only  a  few  representative  simulation  results  are  presented  in  Table  6.  Hie  most  important 
aspect  of  quantile  estimation  in  hydrological  applications  is  estimation  of  the  extreme 
upper  quantiles,  particularly  for  heavy-tailed  GEV  distributions  with  k  <  0.  Table  5 
gives  the  bias  and  standard  deviation  of  the  estimated  upper  quantiles  for  two  GEV 
distributions,  one  with  k  <  0  and  one  with  k  >  0.  Results  are  presented  for  the  ratios 

a  A 

x(F)/x(F)  rather  than  for  the  x(F)  themselves,  since  the  former  quantities  are  more 
easily  compared  at  different  F  values.  For  sample  size  n  >  50  the  three  methods  have 
comparable  performance.  In  small  samples  the  upper  quantiles  obtained  by  MM  estimation 
are  rather  biased,  but  are  still  preferable  to  the  maximum-likelihood  estimators  since 
these  have  very  large  biases  and  standard  deviations.  The  errors  in  the  maximum-likelihood 
quantile  estimators  arise  chiefly  from  a  small  nimber  of  simulated  series  which  yield  large 
negative  estimates  of  k,  and  consequently  give  very  large  estimates  of  extreme  upper 
quantiles. 

Estimation  of  extreme  lower  quantiles  tends  to  be  less  important  in  practice  than 
estimation  of  upper  quantiles,  so  simulation  results  for  this  case  are  not  given  in  detail. 
All  three  methods  give  comparable  results  when  n  >  50,  but  for  small  samples  the  PWM 
estimators  have  smallest  standard  deviation  and  small  or  moderate  bias,  and  are  generally 
to  be  preferred. 
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Table  6.  Biaa  and  standard  deviation  of  estimators  of  GEV  quantiles. 

Tabulated  values  are  bias  and  standard  deviation  of  the  ratio 
* 

x(F)/x(F)  rather  than  of  the  quantile  estimator  x(F)  Itself. 


P-0 

f 

1 

I 

1 

t  o\ 

1  s 

-  It— 0 

P-0. 

1.2 - 

99 

P«0. 

.999 

P-0 

.9 

—  k-0.2  - 

P-0.99 

F-0.999 

x(P)- 

2.84 

x{  P)  - 

■7.55 

x(F)“14.90 

x(F>  — 

1.81 

x(F)-3. 01 

x(F)-3 

.74 

n 

Method 

bias 

s.d. 

bias 

s.d. 

bias 

B*d« 

bias 

S.d. 

bias  s.d. 

hiss  8 

•  d. 

MM 

-.08 

.32 

-.03 

.49 

.12 

.93 

-.04 

.21 

.09 

.31 

.28 

.57 

15 

ML 

.01 

.39 

.64 

2.78 

• 

• 

-.05 

.21 

.17 

1.49 

• 

• 

JS 

.08 

.39 

.00 

.59 

.06 

1.00 

.04 

.22 

.00 

.31 

.05 

.53 

Ml 

-.03 

.27 

.01 

.45 

.13 

.78 

-.02 

.18 

.06 

.24 

.16 

.38 

25 

ML 

.01 

.30 

.17 

.72 

.59 

1.93 

-.04 

.17 

.04 

.24 

.00 

.38 

JS 

.01 

.29 

-.01 

.45 

.05 

.76 

-.01 

.18 

-.01 

.24 

.02 

.37 

MM 

-.02 

.19 

-.05 

.34 

-.03 

.66 

-.02 

.12 

-.02 

.15 

-.01 

.22 

SO 

ML 

-.02 

.19 

-.  02 

.34 

.02 

.67 

-.01 

.12 

.00 

.16 

.01 

.23 

JS 

-.02 

.19 

-.04 

.34 

.00 

.69 

-.02 

.12 

-.01 

.16 

.00 

.22 

tm 

-.01 

.  14 

.00 

.23 

.03 

.38 

-.01 

.09 

.01 

.12 

.03 

.  16 

100 

ML 

.00 

.14 

.02 

.25 

.08 

.44 

-.01 

.09 

-.02 

.11 

-.02 

.15 

JS 

.03 

.15 

.01 

.23 

.01 

.36 

.01 

.09 

.00 

.12 

-.01 

.17 

*  indicates  values  which  varied  widely  between  different  seta  of  1000 
simulations  and  consequently  could  not  be  eatisuited  accurately. 


All  the  methods  of  quantile  estimation  are  very  inaccurate  when  estimating  extreme 
quantiles  in  small  samples  with  k  <  0.  It  is  of  course  to  be  expected  that  a  quantile 
x(F)  cannot  be  estimated  reliably  from  a  sample  of  else  n  if  P  >  1-1/n.  The  implica¬ 
tion  of  this  fact  for  hydrological  practice  is  that  when  estimating  the  upper  quantiles  of 
a  flood  frequency  distribution  for  a  site  with  scanty  data  one  should  seek  to  incorporate 
information  from  other  nearby  sites.  Such  a  regional 1 nation  procedure  can  be  based  on  the 
Ml  estimation  method  for  the  GEV  distribution!  see  Hosking  et  al  (1984). 


Testing  whether  the  shape  parameter  is  sero 
The  type  I  extreme-value  distribution,  or  Quabel  distribution,  is  a  particularly 
simple  special  case  of  the  GEV  distribution,  and  it  is  often  useful  to  test  whether  a  given 
set  of  data  is  generated  by  a  Gumhel  rather  than  a  GEV  distribution.  This  is  equivalent  to 
testing  tdiethar  the  shape  parameter  k  is  zero  in  the  GEV  diatribution.  A  test  of  this 
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hypothesis  nsy  be  based  on  the  PWM  estimator  of  k.  On  the  null  hypothesis  HQ  :  k  -  0 
the  PWM  estimator  k  is  asymptotically  distributed  as  N( 0, 0. 563S/n)  so  the  test  may  be 

*  y, 

performed  by  comparing  the  statistic  Z  •  k(n/0.5635)  t  with  the  critical  values  of  a 
standard  Normal  distribution.  Significant  positive  values  of  Z  imply  rejection  of  Hg 
in  favour  of  the  alternative  k  >  0,  and  significant  negative  values  of  Z  imply 
rejection  in  favour  of  k  <  0. 

The  size  of  the  test  based  on  Z  for  various  sample  sizes  and  its  power  for  sample 
size  50  are  given  in  Tables  7  and  8.  These  results  are  based  on  computer  simulations 
of  1000  samples  for  each  value  of  n  and  k.  The  results  may  be  compared  with  Hosking's 
(1984a)  survey  of  tests  of  this  hypothesis:  the  Z-test  has  power  almost  as  high  as  the 
likelihood-ratio  test  and  for  samples  of  size  25  or  more  its  distribution  on  Hg  is 
accurately  approximated  by  the  standard  Normal  significance  levels.  Since  the  statistic 
Z  is  very  simple  to  compute,  the  Z-test  can  be  strongly  recommended  as  a  convenient  and 
powerful  indicator  of  the  sign  of  the  shape  parameter  of  the  GEV  distribution. 


Table  7.  Empirical  significance  levels  of  the  statistic  Z  for  testing  the 
hypothesis  H  :  k-0  against  one-sided  and  two-sided  alternatives. 


Alternative: 

k<0 

k>0 

k^0 

Nominal  level: 
Sample  size 

10% 

5% 

10% 

5% 

10% 

5% 

15 

10.5 

4.0 

7.3 

3.4 

7.4 

2.5 

25 

9.8 

4.8 

9.4 

4.9 

9.7 

4.4 

50 

11.5 

5.4 

8.4 

4.7 

10.  1 

4.5 

100 

10.5 

4.9 

10.3 

5.4 

10.3 

4.  5 

200 

10.0 

5.7 

9.0 

4.7 

10.4 

5.2 

500 

10.3 

4.9 

9.8 

5.0 

9.9 

5.6 

17- 


Table  8*  Power  of  the  Z-test  of  the  hypothesis  k*0  against  one-sided  and 

two-sided  alternatives.  Saapla  sice  50.  nominal  signficance  level 

S« 


Alternstivei  k<0  k>0  k^O 

k 


-0.5 

.96 

- 

.94 

-0.4 

.91 

- 

.85 

-0.3 

.78 

- 

.71 

-0.2 

.55 

- 

.45 

-0.1 

.27 

- 

.18 

0.0 

.05 

.05 

.05 

0.1 

- 

.17 

.  10 

0.2 

- 

.49 

.36 

0.3 

- 

.83 

.72 

0.4 

- 

.96 

.93 

0.5 

- 

1.00 

.99 

Conclusions 

Estimators  of  parameters  and  quantiles  of  the  GEV  distribution  have  been  derived  using 
the  method  of  probability  weighted  moments.  These  estimators  have  several  advantages  over 
existing  methods  of  estimation.  They  are  fast  and  straightforward  to  cospute  and  always 
yield  feasible  values  for  the  estimated  parameters.  The  biases  of  the  estimators  are 
small,  except  when  estimating  quantilss  in  the  extreme  tails  of  the  GEV  distribution,  and 
decrease  rapidly  aa  the  sample  sise  increases.  The  standard  deviations  of  the  MM 
eatimators  are  comparable  with  those  of  the  maximum-likelihood  estimators  for  moderate 
sample  sices  (n  •  50.  100)  and  often  substantially  less  than  those  of  the  maximum  likeli¬ 
hood  estimators  for  small  a  am  as  (n  -  IS,  25).  MM  estimators  of  GEV  parameters  and 
quantiles  have  asymptotic  Normal  distributions  and  the  large-sample  approximation  to  the 
variance  of  the  estimators  is  adequate  for  sample  sices  of  50  or  more.  Although  MM 
estimators  are  asymptotically  inefficient  compared  to  maximum-likelihood  estimators,  no 
inefficiency  is  detectable  in  saaqiles  of  sice  100  or  less.  The  MM  estimator  of  the 
shape  par*****7  *  of  the  GEV  distribution  nay  be  used  as  the  basis  of  a  test  of  the 

hypothesis  H0  »  k  •  0,  and  this  test  la  simple  to  perform,  powerful  and  accurate. 
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Appendix  1.  Probability  weighted  momenta  for  the  GEV  distribution 
For  the  GFV  distribution  we  have  from  (2)  and  (3) 

6r  “  M0,r,0  "  JolK  *  0,(1  "  ('log  r>k>A)FrdF 

■  Jq{£  +  ad-u^lAls  (r+1)udu  substituting  u  -  -log  F  , 

,  -  .  ...  ,m  -(r+1)u.  ,  ...  ,«*  It  -(r+1)u. 

«  (C  ♦  a/k)  J0  e  du  -  (a/k)  ]g  u  e  du 

-1  -1+V 

-  (C  +  VkHr+l)  -  (aAXrd)  Fd+k)  provided  that  k  >  -1  , 

-  (r+1)-1{5  +  a{l  -  (r+1)'1‘riHkl)/lt]  .  (A1) 


Appendix  2.  Feasibility  of  PHH  estimates  of  the  GFV  parameters 

The  PWM  estimator  k  satisfies  (13),  and  therefore  k  >  -1  provided  that 

(2brb0)/(3b2-b0)  >  j  . 


Now 


and 


2b 


rbo  ’  J  (xi-x)) 


(A2) 


(A3) 


3b2"b0  “  n(n-1 )( n-2)  1>^>k<2xl"Xj"Xk) 


(A4) 


are  both  positive,  so  (A2)  reduces  to  bQ  -  4b^  +  3bj  <  0.  But  we  can  write 

b0-4b1+3b2 .  n(n.,;(n,2l  1>>j>k(-xj+v  » 

thus  bQ  -  4b 1  +  3b2  <  0  almost  surely  and  therefore  k  >  -1  almost  surely.  Results 
(A3)  -  (AS)  above  are  easily  proved  by  induction  on  the  sample  size  n. 

Furthermore,  since 

.  (2b1-bg)k 


(AS ) 


r(1+k)(1-2'k) 


( A6) 


and  2b!  -  bQ  >  0  as  noted  above,  k/(1-2_b)  >  0  for  all  k  and  r(l-Ht)  >  0  because 


k  >  -1,  It  follows  that  we  must  have  a  >  0 


Aggendix  3.  Asymptotic  distribution  of  the  bf 

The  statistic  br  may  ba  written  as  a  linear  combination  of  the  order  statistics  of  a 

random  samplai  we  have 


b 


r 


c<'V 
nj  j 


(A7) 


where  c{**  »  ( j-1) . . . ( J-r )/f (n-1) . . . (n-r)}  and  x,  «  x_  <...<  x  is  the  ordered 
nj  i  2  n 

sample,  As  n  ♦  •  and  j  ♦  "  with  j/n  ♦  ft,  0  <  6  <  1,  cj^'  is  asymptotically  a 
function  of  tha  plotting  position  j/(n-»1)t  in  fact  c„j*  ~  (l/(n+1)}r.  It  is  straight¬ 
forward  to  verify  that  br  satisfies  the  conditions  of  Theorem  1  of  Chernoff  et  al  (1967), 
and  from  that  theorem  it  follows  that  br  is  asymptotically  Normally  distributed  with  mean 
&z  and  variance 

n"’v  -  2n_1  J  J  {F(x)}r{F(y)}r.F(x){1  -  F(y)}dxdy  .  (A8) 

x<y 

A  similar  argument  applies  to  any  linear  combination  of  the  br,  r  •  0,1,2,...,  and  it 
follows  that  the  br  are  asymptotically  Jointly  Normally  diatributed  with  covariance  given 
by 

*xa  m  11«  n  cov(bf,b#)  “  "5  *9r,  +  9>r>  W) 


where 


gr,  -  2  J  J  {r(x)}r{FCy>)“.P(x){l  -  F(y)>dxdy  .  (A10) 

x<y 

To  evaluate  the  qra  for  the  CSV  distribution  we  consider  first  the  case  V  >  0  and 
let 


ao  that 


Ir,  -  2  J  J  {F(x))r{F(y)}*dxdy 
x<y 


(All) 


*re  "  xr+1,s  ”  1r+1,s+1  *  u12) 

1/k 

Substituting  (7)  in  (All)  and  making  the  further  substitution  u  ■  { 1  -  k(x-€)/a>  , 

v  -  { 1  -  kly-O/a)1^*,  we  have 
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_  2  ,«  ,«•  k-1  -ru.  k-1  -«v . 

I  »  2a  I  _  I  u  e  du .  v  e  dv 

rs  J0  Jv 

»  2a2r~k  Jg  vk  'e  SVf(k,rv)dv 

-  a2k"2(r+B)~2kr(1+2k)2F1  [1,2k»  1+k»  s/(r+s>) 


(Gradshteyn  and  Ryrhik,  1980,  pp.  317,  663  )t  hero  r (•,•>  is  the  incomplete  gamna  function 
and  2Ft  is  the  hypergeometric  function.  It  is  convenient  to  transform  the  hypergeometric 
function  in  (A13),  using  results  from  Gradshteyn  and  Ryshik  (1980,  p.  1043): 


^[1,21(1  1+kj  s/(r+s) ]  -  22k{f(1+k)}2/r< 1+2k)  if  r  -  s 


»  {r/(r+s)}  2kG(s/r) 


if  r  >  a  , 

-{•/(ri-sir^GU/s)  +  2r‘ks‘k(r-*-s)2k{ni+k))2/r(H-2k)  if  r  <  s  , 


(A14) 


where  G  denotes  the  hypergeometric  function  G(x)  “  2F 1  ( k ,  2k  i  1+ki  -x)  -  note  that 
G<0)  -  1.  Substituting  back  into  (A12)  and  (A9)  we  obtain  the  following  expressions  for 


the  vrs, 


v  -  a2k-2 (r+1 )”2k ir ( 1+2k)G(r/(r+1 ) )  -  {r(1+k>}2]  , 


(A15) 


v  -  i  a2k"2[(r+2)"21'r(1+2k)G{r/(r+2)}  + 

x / r+ i  * 


(A16) 


(r+1)”k{(r+1)_k  -  2 .  ( r+2 ) }  { T ( 1  +k ) }  2)  , 


v  -  i  a2k“2[(r+s+1)”2kr( 1+2k)G{r/(r+e+1>>  -  (r+s)”kr ( 1+2k)G{ (r+1 )/(r+s) }  + 

r,r+s  2  (AIT) 

2(r+1)  k{(r+s)”k  -  <r+s+1)  k)  {T(  1+k)}  2) ,  s  >  2  . 

When  k  <  0  the  foregoing  argument  is  not  valid  because  the  integral  in  (AH)  does 
not  converge.  However,  the  expressions  <A15)  -  (A17)  are  analytic  functions  of  k  for 
all  it  >  -  and  hence  by  analytic  continuation  expressions  (A15)  -  (A  17)  are  valid 
solutions  of  the  integral  representation  (A9)  -  (A10)  throughout  the  domain  -  j  <  k  <  ®. 
At  the  value  k  *  0,  the  vrg  are  given  by  the  limits  of  (A155  -  ( A 1 7 )  as  k  ♦  0;  these 
limits  are  well-defined. 
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The  results  stated  in  this  Appendix  are  valid  for  arbitrary  positive  integers  r 
and  s,  though  only  the  cases  r,s  -  0,1,2  are  required  for  deriving  the  asyaptotic 
distributions  of  NK  estimators • 


1 


1 
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